Emergent percolation length and localization in random elastic networks 
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We study, theoretically and numerically, a minimal model for phonons in a disordered system 
that shows rich behavior in the localization properties of the phonons as a function of the density, 
frequency and the spatial dimension. We use a percolation analysis to argue for a Debye spectrum at 
low frequencies for dimensions higher than one, and for a localization/delocalization transition. We 
show that in contrast to the behavior in electronic systems, the transition exists for arbitrarily large 
disorder, albeit with an exponentially small critical frequency. The structure of the modes reflects 
a divergent percolation length that arises from the disorder in the springs without being explicitly 
present in the definition of our model. We calculate the speed-of-sound of the delocalized modes 
(phonons) and corroborate it with numerics. We calculate the critical frequency of the localization 
transition at a given density and test the prediction numerically using a recursive Green function 
method. 

PACS numbers: 63.50,+x, 63.20.Pw, 72.20.Ee, 62.30.+d 



Despite decades of experimental and theoretical ef- 
forts, the mechanical properties of disordered solids re- 
main poorly understood. Long-standing challenges in 
materials science such as explaining the low-temperature 
properties of glasses and ultrasound propagation in gran- 
ular media require a robust understanding of the effect 
of disorder on vibrational modes and their localization 
properties. Interest in phonon localization in disordered 
solids has been rekindled by recent experiments on col- 
loidal glasses that mimic several vibrational properties of 
molecular glasses on larger length and time scales [IH3]. 
A significant advantage of colloidal glasses is that their 
vibrational modes can be experimentally reconstructed, 
and they are found to exhibit fascinating localization 
properties. 

Since Anderson's seminal work [5 , the localization of 
electronic wave-functions in disordered metals has been 
the subject of several successful theoretical investigations 
[6]. The study of the vibrational modes of a disordered 
network of masses and springs (see Fig. [I]), which is anal- 
ogous to the electronic problem, received significantly less 
attention. Previous approaches to study this problem 
include perturbative analysis [7] that applies to weakly 
disordered solids, transfer matrix techniques that are re- 
stricted to one dimension [8UTU]. numerical studies [TT1 - 
[T6] and renormalization group approaches [3] [17] [18j [20] . 
One of the successful approaches to modeling vibrations 
in disordered systems has been using the framework of 
Euclidean random matrices [7J [211429] . This problem is 
mathematically related to the localization problem on 
sparse graphs [3QH33] , to the problem of a random walker 
in a disordered environment [3j |34j [35] , and was also used 
to explain aging and slow relaxations in electron glasses 
[M38]. 




M (b) 
FIG. 1. (a) A schematic example of a network of springs and 
masses. All masses are identical with value arbitrarily chosen 
as unity, while the springs are broadly distributed and their 
value is exponential in the Euclidean distance between each 
pair of masses. In the figure, only springs with a strength 
exceeding a threshold K op t are shown, see Eq. (1) of the SM. 
(b) A sketch of the percolating network, yielding a meshwork 
of pathways, and its characteristic length scale I* [4 . We 
will show that for system size L ^> /*, and for dimensions 
d > 1, one can effectively describe the low- lying modes using 
continuum elasticity, leading to a Debye spectrum. 

Many important questions still remain open regarding 
the nature of the localization transition of the vibrational 
modes in disordered systems, its dependence on the di- 
mensionality and the type of the disorder (non-uniform 
masses versus non-uniform springs), and amorphous ma- 
terials such as colloidal glasses still pose many theoretical 
challenges. Most studies rely on numerical solutions, or 
study a particular asymptotic limit, and analytic solu- 
tions for tunable disorder and system dimensionality are 
scarce. In this work, we present a model that is ana- 
lytically tractable, leading to insights into the way the 
localization transition depends on the disorder magni- 
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tude and dimensionality, and elucidating the emergence 
of a disorder-dependent length scale which determines 
various elastic properties of the system. 

Continuum elastic theories of disordered solids typi- 
cally break down below a mesoscopic length scale I* (con- 
trolled by the amount of disorder) intermediate between 
particle diameter and system size [39ti4Tj . In jammed 
packings of grains and random networks, £* exhibits a 
power law divergence whenever the network connectivity 
or packing density are lowered towards the critical val- 
ues necessary to hold these systems rigid. At low densi- 
ties (weak connectivity) the effect of disorder is stronger. 
Plane- wave vibrational modes, which satisy the Deybe 
density of states, exist only below a critical frequency 
uu* ~ I /I* that goes to zero at the critical point. 

In the following, we present a minimal model for a 
random elastic network, which we treat analytically us- 
ing percolation theory. We find that above two dimen- 
sions there is always a localization transition at finite 
frequency, for arbitrarily large disorder. We find an ap- 
proximate analytical expression for the critical frequency, 
which we verify numerically using a recursive Green func- 
tion and finite-size scaling technique. For the delocal- 
ized modes, we calculate, analytically and numerically, 
the speed-of-sound, which is directly related to the low- 
frequency spectrum. The mesoscopic length scale t* 
emerges naturally from the percolation picture, and plays 
a significant role in the calculation of the speed-of-sound 
and the phase diagram. 

The model- Our model of the vibrational properties of 
disordered solids is schematically illustrated in Fig. [T^i, 
which shows point particles (with equal masses, m = 1) 
randomly and uniformly distributed in a d-dimensional 
space, with average nearest-neighbor distance r nn . Be- 
tween every pair of points {i, j} there is a spring with 
spring constant K^. We choose: 

Kij = ^e"^, K u = -mJ^Kij, (1) 

where r^ is the distance between particles i and j, U 
is a characteristic energy scale and £ is the range of the 
exponential interaction. The vibrational frequencies uui 
are related to the eigenvalues A^ of the N x N matrix 
K, as uuf = — A^. This model describes scalar elasticity, 
i.e., it neglects the vectorial nature of the forces [3| f2Tj . 
The randomly chosen points represent the equilibrium 
locations for the masses. In our simple model it is the di- 
mensionless ratio e = £/r nn that controls the strength of 
the disorder - the low density limit (e « 1) corresponds 
to strong disorder. The sum of every row and column 
of K vanishes, due to the negative diagonal elements. 
This 'sum-rule' is a result of momentum conservation in 
the problem, which has important consequences for the 
localization properties and the resulting vibrational spec- 
trum [42] [43] . This sum-rule is associated with a purely 



delocalized mode with uu = 0, corresponding to a uni- 
form translation of all masses. The vibrational problem 
is mathematically equivalent to that of a classical ran- 
dom walker in a disordered landscape [35] [44]; for the 
diffusion problem K describes the Markov process, and 
the probability of being on site j satisfies pj = Kjipi] the 
sum-rule now enforces probability conservation. 

It is not obvious a priori whether the low-frequency 
vibrational modes of this system are localized. It was 
previously shown that the states with uu > uu* are local- 
ized, with uu* — » at low density [3 j. However, since there 
is always a delocalized mode at uu = 0, from continuity 
one might expect that in a large enough system, modes of 
sufficiently low frequency would also be delocalized. We 
find that this is indeed the case above two dimensions. 
Another, perhaps simpler, question regards the structure 
of the density-of-states (DOS) of the eigenvalue distri- 
bution at low frequencies [45]. Will it be the same as 
that of an ordered system, manifesting a Debye spec- 
trum, P(uu) ~ 00 cd (with c the speed-of-sound), or can 
the disorder change it qualitatively? One may think that 
a Debye spectrum and delocalized modes should go hand- 
in-hand, but we shall shortly give a counter-example in a 
one-dimensional system. In dimensions higher than one, 
however, we will show that one always obtains a Debye 
spectrum at sufficiently low frequencies. 

In one dimension and at low densities, the model we 
study reduces to a chain of equal masses connected by 
random springs. For low densities we can neglect matrix 
elements beyond nearest-neighbors, and obtain a power- 
law distribution of the remaining ones [9] [10] . For suffi- 
ciently small disorder (corresponding to large enough e) 
one obtains a Van-Hove singularity of P(uu) = const at 
small eigenvalues, which is Debye-like, while for densi- 
ties lower than the critical density e = 1 a sharper power 
law P(uu) oc l/u a occurs, with < a < 1. In the one- 
dimensional case, all the eigenmodes are localized for any 
density, with the localization length diverging as uu — » 0. 
The fact that the spectrum can still be Debye-like for 
e > 1 shows that one can have localized states and De- 
bye simultaneously. Note, however, that at high densities 
the nearest-neighbor model and the one discussed here do 
not match, since in our case far neighbors will also be im- 
portant. We now proceed to analyze the model in higher 
dimensions, treating separately the cases of high and low 
densities. 

Structure of the eigenmodes at low density.- A key step 
in our approach is to realize that the vibrational proper- 
ties of this spring network can be effectively described 
in terms of an emergent percolation network (shown 
schematically in Fig. [lj}) constructed by retaining only 
springs with > K opt . The optimal spring constant 
K opt will be shown to be proportional to the critical 
spring constant K c , which is defined to be the maximal 
spring constant such that all the springs in the system 
with values greater than K c form an infinite percolat- 
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ing cluster. If we retain only the springs in the infinite 
cluster, the spring network will span the system, but will 
nevertheless have an infinite compressibility (resistivity 
in the analogous electronic problem) , since a macroscopic 
load will have to be held by a single 'bottle-neck' spring. 
If, on the other hand, we retain all springs larger than 
K, with K < K c , we will obtain not a single percolating 
path, but a grid of such paths, with a typical mesh size 
of order /*, see Fig. [TJd. 

The value of K opt can be found by optimizing the 
compressiblity, as shown in the Supplementary Material 
(SM), leading to K opt = K c e~ v ( d - 2 \ This construction is 
made with complete analogy to the electronic system [4] , 
with inverse compressibility playing the role of conduc- 
tivity. The percolation length Z* signals the scale below 
which continuum elasticity breaks-down, and the scale 
above which fluctuations in the compressibility become 
relatively small [4]. As described in the SM, I* ~ l/e v , 
and thus it diverges for low densities. 

We can 'coarse-grain' the sparse percolation network 
of masses and springs, over a length scale of the order 
of the percolation correlation length P, as illustrated in 
Fig. [I] Since K(r) decays exponentially, the effective 
spring constant of each coarse-grained box can be approx- 
imated by the optimal spring K opt with all other springs 
effectively stiff (K > K opt ) or broken (K <C K opt ) g|. 
By coarse- graining over regions of order /*, the result- 
ing distribution of critical springs will be narrow, and we 
have essentially replaced the disordered spring network 
by an approximately ordered one; the fluctuations in the 
values of the spring constants have been suppressed by 
the coarse- graining procedure. Note that the topology of 
this network is ordered by construction, since we coarse- 
grained over boxes. This reasoning suggests that locally, 
the spatial structure of the low-energy eigenmodes should 
coincide with that of the percolating network, but that 
the coarse-grained network would support plane- wave, 
delocalized modes, as would a perfectly ordered system. 
This procedure is only defined in dimensions greater than 
one, where percolation occurs. 

Using the coarse-graining construction we replaced the 
system by an ordered one, and thus we can readily find 
the speed-of-sound of the delocalized modes at low den- 
sities: c = opt /ml* : , where m is the effective mass 
of each coarse-grained oscillator, and the critical spring 
constant K c can be found from continuum percolation 
theory[4 ; namely, K c ~ e~ p ^ e , with P d = 2{r ]c /V d ) 1 ' d , 
where Vd is the volume of the <i-dimensional unit sphere 
and T] c is the standard continuum percolation threshold 
[13SHL giving P 2 = 1-20 and P 3 = 0.87. As previously 
discussed I* ~ ^7, with the critical exponent v equal to 
4/3 in 2D [49]. The fractal dimension of the percolat- 
ing cluster Dp is smaller than the system dimensionality 

and is approximately 1.9 in 2D and 2.5 in 3D [50] . 
Since the size of each coarse-grained box is comparable 
to P d , we expect the total mass of each block to scale as 




-3 -2-10 1 

loge 

FIG. 2. The speed of sound of the low frequency delocalized 
modes as a function of the density parameter e, for a two- 
dimensional system. A low value of e corresponds to a highly 
disordered, sparse network, where the distribution of spring 
constants is very broad. The energy scale U = 1, see Eq. 
and r nn = 1. The straight line shows the exact prediction 
at high densities, below Eq. Q. The other line shows the 
prediction of Eq. ^ of an exponentially suppressed speed 
at low densities, as predicted from percolation theory, taking 
the prefactor as a fitting parameter. The value of Pd depends 
on the system dimensionality, and it is equal to 1.20 in two 
dimensions [35], in good agreement with the numerical results. 

m ~ e d l* Dp ~ e d-vD p _ Combining these results for P, m 
and K opt gives: 

This form is corroborated numerically in 2d in Fig. [2] 
comparing the exact diagonalization of large matrices. 
For each e < 0.2, we determined c by averaging results 
from 2000 realizations with N = 20000 sites with free 
boundary conditions. For each e > 0.2 we used a single 
realization with TV = 10000. For each matrix, the low- 
lying eigenvectors were determined. We discarded the 
uniform mode and Lifshitz modes [5] , where a small num- 
ber of sites are isolated from their neighbors. The small- 
est remaining eigenvalue — Uq gave c = Lujq/tt, where the 
system has linear size L — ~N x l d /e, with £ = 1. For e from 
0.2 to 0.044, the standard deviation in values of c ranged 
from 1% to 30% of the mean, respectively. The value 
shown in Fig. [2] was well determined, with statistical er- 
ror smaller than the symbol size. The SM shows further 
results that the degeneracies of the low-lying eigenvalues 
match what is expected from plane waves in an ordered 
box geometry, see Table SI. Moreover, Fig. 1 in the SM 
shows the average spatial profile of these low frequency 
modes, which is a plane- wave. 

Following the same reasoning as the standard one 
for ordered systems, the DOS which follows from these 
plane- wave- like modes will be a Debye spectrum, and at 
low frequencies we find that: P(oS) ~ ^d~- In the SMs 
we give an alternative derivation for the form of the DOS 
at low frequencies, which relies on mapping this problem 
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to the diffusion problem of a classical random walker in 
a random landscape. Since the prefactor of the Deybe 
DOS depends on the speed-of-sound, this mapping also 
provides an additional way to find the exponential de- 
pendence of c on 1/e, consistent with Eq. ([2|. 

The delocalized plane-wave modes should cease to ex- 
ist as their wavelength approaches the coarse-graining 
scale P, so we can estimate the derealization transition 
frequency 



-2r i i 



CJ* = ck* 



I" 



= e 



vD p /2-d/2 



(3) 



We have used a recursive Green function and finite-size 
scaling technique to systematically investigate whether 
states are localized or delocalized in the infinite-size sys- 
tem [6j [9] , with results shown in Fig. [3j We adapted the 
technique from the one used in Ref. [54] which studied 
the same system without the diagonal terms of K. See 
the SM for details. The percolation prediction of Eq. [3] 
indicates that the delocalized modes should persist to ar- 
bitrarily small e. Due to numerical precision limitations, 
we cannot study e < 0.05 in three dimensions and thus 
cannot numerically determine if there is a smaller e at 
which all states are localized. The prediction of Eq. [3] 
is in good agreement with the numerically determined 
boundary, as shown in Fig. [3j for both two and three 
dimensions. The numerical results for 2d also seem to 
suggest a phase transition, in contrast to the behavior of 
electronic systems, where it is known that arbitrary small 
fluctuations will lead to localization in the orthogonal 
ensemble (which is the appropriate ensemble considering 
the symmetries of our problem). Further numerical and 
analytical study is needed, however, in order to firmly es- 
tablish the existence of delocalized eigenmodes in d = 2. 

High density.- For e ^> 1, the behavior is reminiscent 
of ordered systems: one can show that plane waves are 
approximate solutions to the eigenvalue problem [21] for 
\k\r nn <C 1 and e ^> 1, where k is the wave vector. It is 
found that: 



J drK(f)[e 



ik-f 



!]■ 



(4) 



If we denote the d-dimensional Fourier transform 
of K(r) by K(k), this leads to [21]: A = -uj 2 = 

[K(k) — K(0)]. For K(r) decaying exponentially 

as in Eq. ([!]), in two dimensions, we find for k£ <C 1 
acoustic phonons with a well-defined speed of sound 
c ~ 1% = V37re 2 £7. Fig. [i] shows the good agreement 
of this prediction with the result of exact diagonaliza- 
tion. In 3D, this argument produces acoustic phonons at 
small fc, with c 2 = WirUe 3 . 

Summary.- We have studied the vibrational spectrum 
of a disordered system, using analytical arguments rely- 
ing on percolation theory and a coarse-graining RG-type 
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FIG. 3. Phase diagram showing the localized and delocalized 
modes dependence on density (e) and frequency in two and 
three dimensions. As predicted in Eq.[3] there are both local- 
ized and delocalized modes at any density, and at low densities 
the critical frequency cj* below which modes are delocalized 
depends on e according to Eq. [3] The only fitting parameter 
in each line is the prefactor, which is found to be of order 
unity. White squares indicate parameters where states are so 
strongly delocalized that numerical results did not converge; 
see SM. 



analysis, numerically exact diagonalization, and finite- 
size scaling. We find that for d > 1, there is a critical fre- 
quency uo* vanishing with density, below which a Deybe 
DOS is observed. This frequency marks a localization- 
delocalization transition of the vibrational modes, which 
we study numerically in two and three dimensions. In 
two dimensions further studies are needed in order to 
ascertain whether this is a true phase transition. We 
analytically account for the speed-of-sound of the delo- 
calized modes, both for the low and high density regimes, 
using a percolation approach. A lengthscale P emerges, 
which affects the speed-of-sound and the phase transi- 
tion. In the future, it would be interesting to study the 
heat transport in these systems, and to compare the re- 
sults of the minimal model presented here with those of 
amorphous solids. 
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Supplementary Material 



FINDING THE OPTIMAL SPRING CONSTANT 
K opt AND THE COMPRESSIBILITY 

We denote by K c the value of the critical spring associ- 
ated with percolation, i.e., it is the largest value possible 
for which springs with spring constants larger than K c 
percolate through the infinite system. As described in 
the main text (and illustrated in Fig. 1), we retain all 
springs with spring constant larger than some value K, 
with K < K c . This will lead to a a grid of such paths, 
with a typical mesh size of order P, see Fig. lb of the 
main text. From percolation theory Z* ~ (r — r c )~ v , with 
r and r c the distances associated with the springs K and 
K c . As in the analogous electronic problem (where re- 
sistivity replaces compressibility), this will lead to an in- 
verse compressibility ft -1 ~ K/l* d ~ 2 (this can be seen 
by replacing the system by an ordered grid with lattice 
spacing I* and the spring constant of each edge in the 
lattice is K). In our problem r = £ log (If), and thus for 
a given choice of K the associated compressiblity is: 

k- 1 ~K[t\og{K/K c )Y. (1) 

As expected, the inverse compressibility vanishes at 
K C1 and thus clearly there is a maximal compressibility 
for some K opt < K c . A straightforward calculation gives: 

K opt = K c e-^ d - 2 \ (2) 

Since K opt oc K c , the scaling of I* is still I* ~ \je v . For 
low densities (small e), corresponding to high disorder, 
this length scale diverges. 

This leads to an inverse compressibility of: 

k~ x (xK c e^ d ~ 2 \ (3) 

The polynomial correction term e v ^ d ~ 2 ^ of the inverse 
compressibility is consistent with the renormalization 
group results on a related model [TJ |2], which suggest 
that the results are robust to the precise details of the 
disorder used. 



DENSITY-OF-STATES AT LOW FREQUENCIES 

We argue that the density-of-states will correspond to 
a Debye spectrum at low enough frequencies and for 
d > 1, for arbitrarily small e. To show this, it will 



be useful to employ the mapping between the vibra- 
tional problem studied here and that of the classical ran- 
dom walker, discussed in the main text. For the diffu- 
sion problem, the return probability Ri(t) is defined as 
the probability to be at site i at time t giving that we 
were at site i at time t = 0. It can be shown that [3 
R(t) = Ri(t) = j P(\)e- Xt d\ where R(t) is the average 
over all sites z, and P(X) is the density of states of the ma- 
trix K. This formula implies that a higher DOS at small 
frequencies corresponds to a larger return probability at 
long times, which suggests slower diffusion. In general, 
a Debye spectrum (in any dimension) corresponds to the 
case of normal diffusion and the return probability scales 
as \jDt d l 2 . Thus, to establish the existence of a Debye 
spectrum at low frequencies, it suffices to prove that at 
asymptotically long times a particle will undergo normal 
diffusion. 

To perform this step, we employ another mapping, be- 
tween random walks and electrical networks. Consider 
a random resistor network whose resistance between two 
sites is proportional to the matrix elements of if^-, i.e., 
the hopping rate. Then, the Einstein relation tells us 
that the diffusion coefficient of the network (if it is non- 
zero) will be proportional to the conductivity. Therefore, 
the question of establishing normal diffusion is equivalent 
to proving that a finite conductivity exists in the net- 
work. This result was proved using percolation methods 
[4] , in an identical system where the matrix elements (re- 
sistors) between two sites depend on both their energy 
difference and distance. For the case of very high (infi- 
nite) temperature, the matrix elements depend only on 
the distance, exponentially, as in the case we discuss here. 
Thus, for dimensions higher than one, the result of Ref.[4| 
shows that there is conductivity, and therefore diffusion, 
at asymptotic times. Moreover, this allows us to derive 
the diffusion coefficient, since percolation theory gives the 
dependence of the conductivity on e: a oc D oc e~ p ^ e [J. 

We now explain how this relates to the results of Ref. 
[3], where a strong peak was shown to occur in the DOS 
at very small frequencies. There, it was shown that for 
low densities the DOS can be approximated by: 

= dV d e- V ^^/^ d [\o g {^/2)) d -\ 

UJ 

where Vd is the volume of the <i-dimensional sphere. 

At low frequencies Eq. Q does not conform to the De- 
bye spectrum, giving a more-rapid decay. Thus Eq. Q 
must fail near some frequency cj*, which should approx- 
imately occur at the derealization transition. Eq. Q 
was derived by approximately calculating all moments 
of the eigenvalue distribution to an accuracy of 0(e m ) 
with m > 1. The critical frequency uj* of Eq. (3) in the 
main text is small enough such that the corrections due 
to the Debye spectrum at low frequencies do not signifi- 
cant change any of the moments. Moreover, we can use 
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the form of a;* to calculate the integral of the density-of- 
states described by Eq. Q, from cj* to the upper cutoff 
uo 2 = 2, corresponding to a pair of masses [3]. The result 
it that there is always a finite fraction of localized modes. 
Taking the Debye spectrum below cj*, P(oS) ~ u cd , and 
using the value of c from Eq. (2) of the main text, we can 
estimate the number of delocalized modes, and find that 
it scales as e vd . As one goes to higher disorder (lower 
e) the fraction of delocalized modes decreases, becoming 
negligible at large disorder. 



NUMERICAL VERIFICATION OF 
DELOCALIZED MODES 

We have tested the analytic arguments numerically, 
working with large, sparse matrices, obtained by truncat- 
ing the full matrices: all elements smaller than a cutoff 
value are set to zero. For a small enough cutoff value the 
results are independent of the value of the cutoff, as they 
should be. 

Some of the low-lying eigenmodes consist of Lifshitz- 
like, rare isolated points [5]. Their participation ratio is 
close to unity, and thus they are essentially localized on 
a single point, with a vanishingly small frequency, ex- 
ponentially small in the distance separating this point 
from others. However, the contribution of these events 
to the density-of-states is negligible for a large enough 
system and for small enough frequencies, since it goes 
to zero faster than any power-law (note that this relies 
on the exponential decay of the matrix elements, and 
for other forms decaying faster with the distance this 
assertion might be false). Aside from the rare Lifshitz 
singularities, the low frequency modes are found to be 
delocalized in these finite systems in two and three di- 
mensions. In two dimensions, the lowest two non-trivial 
frequencies are nearly degenerate, and in three dimen- 
sions triply degenerate, as expected from the analysis of 
eigenmodes in an ordered-lattice geometry. Looking at 
higher non-trivial frequencies, their values are consistent 
with the coarse-grained ordered picture, predicting par- 
ticular ratios of the low lying frequencies, as shown in 
Table 1. 

The average spatial profile of these low frequency 
modes clearly shows plane wave oscillations, as shown 
in Fig. [I] 



Index 


lj, 2D ord. 


cj, 2D num. 


cj, 3d ord. 


cj, 3d num. 


1 


1 


1 


1 


1 


2 


1 


1.022 


1 


1.014 


3 


V2 


1.415 


1 


1.033 


4 


2 


1.994 


V2 


1.390 


5 


2 


2.004 


V2 


1.397 


6 


75 


2.222 


V2 


1.426 


7 


75 


2.244 


Vs 


1.678 


8 


Vs 


2.826 


2 


1.950 


9 


3 


3.001 


2 


1.991 


10 


3 


3.042 


2 


2.037 



TABLE I. Comparison of the frequencies of the low lying 
modes of the perfectly ordered case (ord.) and the numerical 
results (num.) on the disordered case in 2d and 3d, for open 
boundary conditions. The first non-trivial frequency is nor- 
malized to 1. For the two-dimensional case N=500,000, and 
for the three-dimensional case for N=100,000, with e — 0.1. 



by stitching together many (d— l)-dimensional strips. In 
the original Anderson model, each site is coupled only to 
its nearest neighbors, so each new slice has direct cou- 
plings only to the closest neighbors; this property is es- 
sential for the recursion. 

In Ref. 7, we extended this technique to non-lattice 
problems, as studied in this paper. Without a lattice, the 
quasi-one-dimensional system is stitched together from 
slices of width w and length /, which contain a random 
number of sites, chosen from a Poisson distribution with 
a mean density e d . In order to have particles in each slice 
interact only with particles in the nearest neighbor slices, 
we must impose a cutoff on the interaction between sites. 



O Numeric 
Cosine prolila 



DETAILS OF THE RECURSIVE GREEN 
FUNCTION CALCULATIONS. 

The recursive Green function technique with finite size 
scaling was first used to study Anderson localization 
by MacKinnon and Kramer [6 . It studies quasi-one- 
dimensional systems with one very long axis and d — 1 
shorter axes, of width w. The system is built recursively 



FIG. 1. The lowest frequency mode is shown for a single ma- 
trix with N = 500, 000, in two dimensions, with e = 0.1. The 
amplitudes were binned along one of the sides of the square in 
which the points are randomly chosen, and the graph shows 
the sum of amplitudes in each bin. While the eigenmode itself 
is not a simple plane wave, after 'coarse-graining' in this way 
one obtain a cosine profile, precisely as would be obtained in 
an ordered system. 
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- co=3.16e-05 / 8=0.097 

- 03=5.626-05, 8=0.097 
-03=1.006-04, e=0.097 

co=1.78e-04, 8=0.097 

- co=3.16e-04, 8=0.097 
co=5.62e-04, 8=0.097 

-co=1.00e-03, 8=0.097 
-co=1.78e-03, 8=0.097 



FIG. 2. Part of the data used to produce the phase diagram 
in Figure 2 of the main text. For e = 0.097 in d = 2, we show 
A w (w). The phase boundary between the localized (sloping 
down) and delocalized (sloping up) energies is clear. 



That is, we set exp(— r) — » for r > L c for some cutoff 
L c . Since we use periodic boundary conditions in the 
d — 1 dimensions of the slice, we choose the cutoff to be 
L c = min[Z, w/2]. Since there is no lattice, w and I must 
be chosen appropriately large so that the system is not 
disconnected. 

We begin with a single slice of width w and length I. 
We add a second slice and find the portion of the energy- 
resolved Green function G(A) that connects any site in 
slice 1 to any site in slice 2. We recursively add more 
sites, always finding the portion of G connecting sites 
in slice 1 to sites in slice TV, G? N (\) = (1\G N (X)\N). 
G± N (X) is an Mi x Mn matrix, where Mi is the number 
of sites in slice i. The probability for a particle injected 
with energy A into slice 1 to make it to any site in slice 
N is 



N 



Tr[G? N W} 1/2 - 



(5) 



Since all states are localized in Id systems, P/v decays like 
exp(—2Nl/£ w ), for localization length £ w . We extract the 
effective localization length at energy A and width w by 



lim 



log -P/v 
Nl ' 



(6) 



N must be chosen large enough to get good statistical 
accuracy for £ w [3|9]. Statistically, the system behaves as 
though it consists of Nl/£ w independent and identically 
distributed samples chosen from a normal distribution, 
so the sampling error decreases as TV -1 / 2 . 

When £ w has been acquired for a range of widths w, the 
theory of single-parameter scaling is used to extrapolate 
to the infinite-size system [6j |9j [10] . According to single- 
parameter scaling, for sufficiently large w that irrelevant 



corrections can be neglected, the dimensionless scaling 
variable A w = £ w /w obeys a universal scaling law 



A u 



0(6, A) 



(7) 



for some unknown universal function F, where 0(e, A) 
characterizes the infinite-size system: it is the local- 
ization length for localized systems and the correlation 
length for delocalized systems. In particular, if A w in- 
creases with w, then the state is delocalized and if A w 
decreases with w, then the state is localized. This is the 
approach adopted to determine the phase diagram in Fig- 
ure 3 of the main text. Single-parameter scaling has been 
demonstrated to hold in similar systems without the sum 
rule [7] , but it has not yet been proved in this system. It 
is possible that the emergence of the percolation length I* 
would cause a breakdown in the single-parameter scaling, 
but the data so far do not show any such signs. 

The details of setting up the recursive Green func- 
tion calculation have been extensively given elsewhere 
[TTj H2] . The basic principle is that, in each step, the 
exact Green function of the wire with N slices G N and 
the exact Green function of the new slice g are known. 
We write Go = G N + g. The perturbation that connects 
them V can be treated non-perturbatively using a Dyson 
equation 

G N+1 = G + G N+1 VG = G + G VG N+1 . (8) 

Matrix elements taken on Eq. [8] give the recursion rela- 
tions required. 

For our case, the sum rule in the definition of the ma- 
trix K means that the perturbation V contains not only 
off-diagonal terms connecting sites in the two slices, but 
also contains diagonal terms, adjusting the on-site en- 
ergy of each site. This does not pose any difficulty for 
the method, but requires more matrix operations per cal- 
culation, increasing the cost of the calculation. 

In each calculation, we choose the values of A we wish 
to study, and for each A store the matrices G± N , Gf^ and 
G^n- For the N + 1st slice being added to the system, it 
has a Hamiltonian K and the off-diagonal perturbation 
connecting the existing wire to the new slice is U . Then 
in the usual problem, without the sum rule, the result is 



^1,1 



(X — K 

Gl,AT 



u^G N 



N.N 



uy 



UKjf N+l,N+n 



r N TT(r N+1 ^ 

Lr 1N U{Lr 1N ) , 



(9) 
(10) 

(11) 



where we used g~ x = A — K '. When we include the sum 
rule, the perturbation contains the off-diagonal terms U, 
the diagonal terms Un acting on the sites in the Nth 
slice, and the diagonal terms Un+i acting on the sites of 
the new slice. Then we find 
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r N+l 
^N+l 



iV+l,iV+ 



n N 



-l n N 



1 = [X-K-U^t MN -G%U N , 

t „ + Un{^-m n — Gn,nUn) 1 G% iN ]UGnXi i n+ii 

Gi,? 1 = G i t i + |G^ Ar C/(G^ 1 ) t + Gi N U n (1 Mn - G^ >Ar C/Ar)" 1 [(G^ Ar ) t + G^^UiG^ 1 ^]^ . 



(12) 
(13) 
(14) 



If the calculations are done all at once, there is no 
need to store Gi,i. But it is often convenient to break 
the calculation of a very long wire into many smaller 
pieces, and then stitch them together. This stitching 
requires that Gi ? i be stored. Consider one wire with N± 



slices with recursively calculated Green functions G 1 \, 
G\n 1 ? an d N and a second wire with N 2 slices with 
recursively calculated Green functions g^J, gi 2 N , anc ^ 
N . Let N = Ni + Then we can stitch g onto G 
to find 



G 



N 

N.N 



„N 2 



iV 



M A r 1 



n N 

{ ^N 1 + 1,N 2 

D 



T N 1 ,N 1 



J Nl+1 g^)-^G^ Ni 



JVi+ifli,i 

[i Mwi -^(K + .- p ' 

Gft = Gft + G± 1 Ni D [ug»X - EWsft]" 1 ^ + ^ } G^ . 



i/l,JV 2 ' 



(15) 
(16) 
(17) 
(18) 
(19) 



It is generally desirable to let N be large enough that 
the statistical error in i w is less than 1%. For values of A, 
e where the states are highly delocalized, i w can be very 
long, requiring an inordinate length N to get reasonable 
errors. Conveniently, the values of £ w are reasonable for 
the localized states and the nearby delocalized points, so 
we can extract the derealization boundary cj*(e), while 
we cannot quite prove that the strongly delocalized states 
are actually delocalized. Care must be taken to avoid 
both overflow and underflow in storing Gi ? ]\r 5 as detailed 
in Ref. [12]. 

Figure [2] shows some of the calculations used to create 
Figure 2 of the main text. The statistical error bars are 
smaller than the data points for all but the lowest energy 
(most delocalized) state. For these simulations, I = 50 
and N ranged from 2 • 10 5 to 2.6 • 10 T . 
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